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With an underlying common theme of competing length scales, we study the many-body Schrodinger equa- 
tion in a quasiperiodic potential and discuss its connection with the Kolmogorov-Arnold-Moser (KAM) problem 
of classical mechanics. We propose a possible visualization of such connection in experimentally accessible 
many-body observables. Those observables are useful probes for the three characteristic phases of the problem: 
the metallic, Anderson and band insulator phases. In addition, they exhibit fingerprints of non-linear phenomena 
such as Arnold tongues, bifurcations and devil's staircases. Our numerical treatment is complemented with a 
perturbative analysis which provides insight on the underlying physics. The perturbation theory approach is 
particularly useful in illuminating the distinction between the Anderson insulator and the band insulator phases 
in terms of paired sets of dimerized states. 
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I. INTRODUCTION 

Ultracold atoms are emerging as a versatile arena for the 
study of a variety of problems in physics. These clean and 
highly controllable systems can be viewed as simulators of 
complex quantum phenomena with applications in condensed 
matter, quantum optics, atomic and molecular physics and 
nonlinear dynamics, as well as in particle physics and cos- 
mology 1 1 , 2 |. Some examples of this trend are the laboratory 
realization of the superfluid to Mott insulator transition with 
bosons 1 3 1, observation of the corresponding metal to Mott in- 
sulator transition with fermions |4, 5 | and the creation of the 
Tonks-Girardeau gas in one dimensional systems |(6lEl. Ul- 
tracold atoms have also provided a laboratory realization of 
Wannier- Stark ladders HI and the kicked-rotor model 121, a 
paradigm in the study of classical and quantum chaos ifTOl . 
Recent successful loading of atoms in a bichromatic optical 
lattice geometries 1 11 1, has opened new avenues for the study 
of phenomena where the competition between various length 
scales is at the heart of the problem. 

Systems with competing lengths have fascinated physicists 
as well as mathematicians in view of their exotic fractal char- 
acteristics 1 12, 13 1 . Such systems with two competing period- 
icities, commonly known as almost-periodic or quasiperiodic, 
occur very commonly in nature. The most commonly studied 
example in quantum physics is the single particle Schrodinger 
equation in the presence of a quasiperiodic (QP) potential. 

Motivated by experimental realization of two-color lat- 
tices, we revisit the problem of the Schrodinger equation 
in a quasiperiodic potential and its relationship with the 
Kolmogorov-Arnold-Moser (KAM) 1 14] problem of classical 
mechanics. Our particular focus is the treatment of many- 
body effects that are present in QP systems that have been 
extensively studied at the single-particle level. By focusing 
on experimentally accessible observables such as the momen- 
tum distribution and density-density correlations, we demon- 
strate the possibility of experimental visualization of the rela- 
tionship between the metal-insulator transition and the KAM- 
Cantori transition 1 15] in many body systems. These observ- 



ables are found to exhibit fingerprints of various paradigms of 
nonlinear systems such as Arnold tongues, bifurcations and 
devil's staircases. One of the aims of this paper is to com- 
municate the excitement of ultracold atomic physics to the 
nonlinear dynamics community. We also hope that cold atom 
community will also benefit from our discussion of the rela- 
tionship between problems of condensed matter theory and 
nonlinear dynamics. 

The many-body systems that we treat here are ensem- 
bles of ultracold spin-polarized fermionic atoms confined in 
one-dimensional bichromatic optical lattices. We treat cases 
in which the two lattices have incommensurate periodici- 
ties, and thus constitute a quasi-periodic potential for atomic 
motion |16|. Since the metal-insulator transition in single- 
particle QP systems is associated with a localization transi- 
tion of extended single-particle states, a natural probe of such 
transitions in many-particle systems are the density distribu- 
tion and density-density correlations. At present, the den- 
sity distribution of atoms confined in a lattice has not been 
easily accessible to experimental measurement, instead the 
quasi-momentum distribution and its corresponding correla- 
tion functions, have been measured after the lattice potential 
has been suddenly removed. 

Using the momentum-position duality of our basic model 
ifTTl [T8 L we show that time of flight images encode local 
density information and can be used to identify the possible 
phases of ultracold atomic gases, e.g., metallic phases and 
various types of insulators. Many of our results are explained 
analytically using perturbation theory, although the complete 
picture is obtained by exact numerical calculations, particu- 
larly near the metal- Anderson insulator transition. 

The paper is organized as follows. Section II contains an 
overview of experiments on ultracold atoms in optical lat- 
tices. Section III describes the basic Hamiltonian, and the 
experimental observables. In section IV, we describe the rela- 
tionship between the metal- Anderson insulator transition and 
KAM-Cantori transition. In Section V we discuss possible ex- 
perimental manifestations of effectively nonlinear behavior in 
these many -body systems. Sec. VI provides a summary of our 
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results and states our conclusions. 



II. ULTRACOLD ATOMS IN OPTICAL LATTICES 



An optical lattice is created by the interference of counter- 
propagating laser beams that give rise to a spatially periodic 
intensity pattern. The intensity pattern corresponds directly 
to a potential for neutral atoms via the a.c. Stark shift of the 
atomic energy levels. The two important parameters of an 
optical lattice are the depth of the lattice potential wells and 
the lattice constant, a. The well depth of the lattice can be 
tuned by changing the intensity of the laser, while a can be 
tuned by changing the wavelength of the laser or by changing 
the relative angle between the two laser beams 1 1 1. 

The systems of interest here are gases of ultra-cold spin- 
polarized fermionic atoms trapped in the lowest band of a 
transversal 2D optical lattice. For most such systems of 
current experimental interest, where atomic interactions are 
short-range, spin-polarized fermions are effectively noninter- 
acting due to Pauli blocking. The 2D lattice depth is made 
strong enough to freeze the motion of the atoms transversally, 
creating an array of independent ID tubes |[6j[T9l. Along the 
axis of the ID tubes, additional optical lattices can be im- 
posed, as has been done in Refs. l7l|2Q|. We discuss cases 
in which two such lattices are imposed, with incommensurate 
periods O (when one of these lattices is much stronger than 
the other, we refer to it as the primary lattice). The combined 
lattices therefore generate an effective quasi-periodic (QP) po- 
tential along the axial direction. Such an experiment has been 
recently implemented for bosonic atoms in Ref. |11|. Our 
treatment considers cases where the ratio of the two lattice 
constants is equal to the "golden mean", a = (>/5 — l)/2, 
which is one of the best-studied examples in single-particle 
physics ifTSl . 

In most experiments atoms are first trapped and cooled to 
quantum degeneracy. They are subsequently loaded into the 
lattice by adiabatically turning on the lattice laser beams. At 
the end of each experimental sequence atoms are probed by 
using time of flight images. These are obtained after releasing 
the atoms by turning off all the confinement potentials. The 
atomic cloud expands and then photographed after it enters 
the ballistic regime. Assuming that the atoms are noninteract- 
ing from the time of release, properties of the initial state can 
be inferred from the spatial images |21|: the column density 
distribution image reflects the initial quasi-momentum distri- 
bution, and the density fluctuations, namely the noise correla- 
tions, reflect the quasi-momentum fluctuations. These quan- 
tities, which will be defined below -see Eqs. ( |5|6| )-, have 
been shown to be successful diagnostic tools for character- 
izing quantum phases and have been recently measured in 
bosonic quasi-periodic systems I 22J . 



III. OUR MODEL SYSTEM: THE HARPER EQUATION, 
MANY-BODY OBSERVABLES AND SELF-DUALITY 

A. The Harper Equation 

If the intensity of the secondary lattice is much weaker 
than that of the primary lattice, the low energy physics of the 
fermionic system can be well described by the tight-binding 
Hamiltonian [2J: 

H = -J ^{c]cj^i + 4+1^^-) + E ^J^J^ 

j 3 

where Cj is the fermionic annihilation operator at the lat- 
tice site j, and J is the hopping energy between adjacent 
sites. The main effect of the QP potential is to modulate 
the on site potential. It is accounted for by the term Vj = 
2Vq cos(27r(Tj + (/)). The parameter Vq is proportional to the 
intensity of the lasers used to create the secondary lattice (161, 
(J is the ratio between the wave vectors of the two lattices 
which we choose to take value a = (>/5 — l)/2, and (j) is 
a phase factor that is determined by the absolute registration 
of the two lattices. 

To model real experimental conditions, averaging over (j) 
is required. This averaging takes into account, on one hand, 
the phase fluctuations from one preparation to another. Those 
arise due to the difficulty to lock the position of the cloud over 
several shots. On the other, the fact that typical experimental 
set-ups generally consist of an assembly of one-dimensional 
tubes with different lengths and phases among them. 

For a single atom, the eigenfunctions V^^-^^ and eigenener- 
gies E^^^ of the Hamiltonian in Eq.([T]) satisfy: 

(2) 

Where A = Vq/J, e^^) = £;(^)/J, and e^^) < e^^+i). Eq. 
^ is known as the Harper equation, a paradigm in the study 
of ID quasiperiodic (QP) systems [18J . For irrational a, the 
model exhibits a transition from extended to localized states 
at Ac = 1. Below criticality, all the states are extended Bloch- 
like states characteristic of a periodic potential. Above criti- 
cality the Harper model becomes equivalent to a correspond- 
ing Anderson model, the spectrum is a pure point spectrum 
and all states are exponentially localized. At criticality the 
spectrum is a Cantor set and the gaps form a devil's staircase 
of measure unity 1 12 |. 

In our numerical studies, a is approximated by the ratio of 
two Fibonacci numbers Fm-i/Fm, (Fi = Fq = = 
Fi -\- Fi-i), which describe the best rational approximant by 
continued fraction expansion of a. For this rational approx- 
imation the unit cell has length Fm and the single-particle 
spectrum consists of Fm bands and Fm — 1 gaps. The gaps 
occur at Qn/2, ±(7r - |(5n|/2) with Qn = ±{27r / a) (ncr) 
reciprocal lattice vectors constrained in the interval aQn G 
(— TT, tt]. Here (na) = na (mod 1), n an integer . We denote 
by Np the number of atoms, Ni = Fm is the number of lattice 
sites, and the filling factor u = Np/Ni ranges from to 1. 
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B. Many-body Observables 

An ensemble of spin-polarized fermions at zero tempera- 
ture are "stacked up" into the single-particle eigenstates of 
increasing energy, with one particle per quantum state. The 
energy of the highest occupied level, which depends on the 
filling factor is the Fermi energy, Ep. Since at the critical 
point all the single-particle wave functions become localized, 
at the many-body level polarized fermions also exhibit a tran- 
sition from metal to insulator at Ac. However, in addition to 
these two phases the fragmentation of the single-particle spec- 
trum in a series of bands and gaps introduces additional band 
insulator phases when the Fermi energy lies within a gap. The 
most relevant insulating phases occur at the irrational filling 
factors: v = Fm-i/Fm and v = Fm-2/Fm ( which respec- 
tively correspond to = a, 1 — a), associated with the leading 
band gaps. In the band insulator phases the many body sys- 
tem is an insulator, irrespective of the value of A, even though 
extended single particle states are occupied. 

Since the metal-insulator transition is clearly signaled by 
the onset of localization of extended single-particle states, a 
natural probe of this transition is the many-body density pro- 
file, pj, and the density-density correlations, A(ji, j2), which 
can be written in terms of single particle wave functions as: 



Pj 

A(jl,j2) 



m—l 



(m) 



El (m)* / (m) 



(3) 



(4) 



Here fij = CjCj and we have used Wick's theorem to evalu- 
ate A(ji, j^). However, in general, such local observables are 
hard to measure experimentally due to the lack of address- 
ability of individual lattice sites for typical lattice spacing. 
Instead, time of flight images access non-local observables 
such as the quasi-momentum distribution, n{Q) and the quasi- 
momentum fluctuations, A((5, Q^), which are given by: 



A((5) = A((5, 0) can be written as: 



rik 



El (m)* (m) 
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A{Q{k)) = E ho 



(m) 



m=l 



E(m)* (m) 



(8) 



. (9) 



C. Self-duality 

The self-duality of the Harper equation corresponds 
to the property that single particle eigenstates ipj and 
their corresponding Fourier transformed eigenstates, Zn = 
_^^-in(f) ^.^-ij{27ran-\-e) satisfy the same equation with 

the roles of J and Vq interchanged ifTTlfTSl : 



2 e 

{Zn-l + Zn^l) + - COs(27rcrn + 0)Zn = -J^n- 



(10) 



This relation implies that the experimentally relevant vari- 
ables, r]k, (see Eq. ^) also satisfy 



-(^Mod[fe+FM-i,FM]6 

2 

-- cos{27rk /FM)r]k - 



■ VModik-FM-i^FM]^^"^) 



-rjk 



(11) 



As discussed below, self-duality is a key to obtaining exper- 
imental information on local quantities from measurements. 
In this context, it is important to understand the relationship 
between the index j of ipj that satisfies the Harper equation 
^ and its corresponding index k in r]k that satisfies the dual 
equation ( pTj ). An example of the mapping is provided in Ap- 
pendix A. It should be noted that Fibonacci sites in real space 
are mapped to Fibonacci sites in the momentum space, up to 
a common displacement. This shift is dependent on the phase 
factor (j). We obtain this relationship numerically by diagonal- 
izing the Harper equation for a given A, labeling the states in 
increasing order in energy, and then finding the corresponding 
momentum space dual by repeating the same procedure but 
with A replaced by 1/ A. 



^(Q) = ]^Ee"^"^'"'^49, (5) 
A{Q,Q') = {h{Q)h{Q'))-{h{Q)){h{Q')), (6) 

where Q is the quasi momentum which can be expressed in 
terms of the indices A: = 0, 1, 2, . . . as: Q{k) = ^ A' ^^^^ 
is constrained to the interval Qa G (— tt, tt] . 
Introducing r]k, the Fourier transform of ipk, 

the observables Uk (also denoted sometimes as n(Q(/c))) and 



IV. LOCALIZATION TRANSITION AS A KAM-CANTORI 
TRANSITION 

The perturbative treatment of the QP potential is fundamen- 
tally related to the treatment of a non-integrable perturbation 
applied to an integrable Hamiltonian system in classical me- 
chanics. Both cases exhibit the small-divisors problem related 
to the presence of high-order terms with small denominators 
in the perturbation expansion. In the classical system, Kol- 
mogorov, Arnold and Moser (KAM) |14| solved the prob- 
lem and demonstrated that most of the invariant tori in the 
phase space are not destroyed by a sufficiently weak nonin- 
tegrable perturbation. Outside the perturbative regime, KAM 
tori break, becoming an invariant cantor set, known as Cantori 
L15J. 
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The study of the Harper equation is fundamentally related 
to the KAM type problems of classical mechanics. Both sys- 
tems share the mathematical difficulty of having higher-order 
terms with small denominators when the quasiperiodic poten- 
tial or non-integrable term are treated perturbatively. Under 
this point of view, the metallic phase in the Harper Equation, 
with continuous spectrum and B loch-type wave functions has 
been identified as the analog of the KAM phase with invariant 
tori, while the localized phase of the Harper system with point 
spectrum and exponentially localized states has been usually 
compared with the Cantori phase. At the single-particle level 
this connection has been visualized by means of the so called 
Hull function ||23l[24]| defined as ipn = e'^^^xi'^^)^ with a a 
real phase factor. x(n(j) is a smooth and continuous function 
in the extended phase but becomes discontinuous for A > 1. 
Here we propose instead to look at the return map of the local 
density of the atomic cloud, pj vs pj+i, as a cleaner visualiza- 
tion of the KAM to Cantori transition at the many-body level. 
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FIG. 1: (color online) Return maps in real space for (a) and (c) 
A = 10 (red), 1 (green), 0.5 (blue) and the corresponding return map 
in momentum space (b) and (d) A = 0. 1 (red), 1 (green), 2 (blue). 
The filling factors are = 0.25 for the upper panels and u = 0.5 for 
the lower panels. The insets are averaged over 50 random phases. 

Fig[T] shows the return map for various rational filling fac- 
tors v = Np/Ni = 1/4, 1/2 and different values of A. For 
A < 1, the return maps are smooth curves and correspond to 
the KAM tori of the extended or metallic phase. For A > 1, 
the density profile is a discontinuous function, a Cantori. The 
discreteness of the return map for generic filling factors can 
be easily understood in the A ^ oo limit, where the wave 
functions are fully localized and thus the return map can only 



take the four possible values (0, 0), (0, 1), (1, 0), (1, 1). Ex- 
actly at the transition point A = 1, the smooth curves become 
disconnected. 

Using the duality transformation, similar return maps can 
be drawn in momentum space, vs nk±ko (^o = ^m-i)- 
While nearest neighbor sites are connected in position space 
due to a finite J, quasi-momentum components separated by 
the main reciprocal lattice vectors of the secondary lattice are 
connected by a finite Vq. We find the momentum return maps 
exhibit an important advantage compared with the local den- 
sity maps, which is related to the fact that they are insensitive 
to variations of the phase and retain their pattern when aver- 
aged over it. This is not the case in the density maps as shown 
in the insets of Fig. [T] 

We now consider the case where the filling factor is a ratio 
of two main Fibonacci numbers (u = a, 1 — a). We refer 
to those filling factors as irrational filling factors since they 
approach an irrational number in the thermodynamical limit. 
As discussed earlier, this results in a band insulating phase as 
the Fermi energy Ep lies in the gap, just outside a filled band. 
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FIG. 2: (color online) Return maps for the irrational filling factor jy — 
a. Panel (a) is in real space with A = lO(red), 1 (green), 0.5(blue) 
and panel (b) is in momentum space with the dual values A = 
0.1 (red), 1 (green), 2 (blue). The insets are averaged over 50 random 
phases. 

In contrast to the generic or rational filling factors discussed 
earlier, the return map for the irrational filling case is found 
to remains smooth regardless of the value of A (See Fig|2]). 
This result may appear somewhat counterintuitive, because as 
A ^ oo, all the single particle wave functions become local- 
ized. However, an exception to this simple picture occurs near 
irrational filling, when the density function exhibits a continu- 
ous distribution. This difference between the rational and the 
irrational filling is due to the presence of a group of paired 
states centered around the most dominant band edges ( asso- 
ciated with the dominant gaps). Each of those states is dimer- 
ized, by which we mean that they localize at two neighbor- 
ing sites. The existence of a pair of such dimerized states is 
the key to understanding the difference between rational and 
the irrational filling: although a single dimerized state causes 
delocalization near band filling (explaining irrational filling 
case), the existence of a pair can cancel the delocalization ef- 
fect. Technical details of this argument are presented in Ap- 
pendix B, where we show that the smooth character of the 
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map at the irrational fillings, a, cr^ associated with the leading 
gaps can be understood by the breakdown of non-degenerate 
perturbation theory in the large A limit. 

In summary, the band insulating phase belongs to the KAM 
phase irrespective of the strength of the disorder. 



V. FINGERPRINTS OF NON-LINEAR PHENOMENA IN 
MANY-BODY OBSERVABLES 

We now study the quasi-momentum distribution of polar- 
ized fermions, which is directly accessible in the time of flight 
images of ultracold atoms. In the following we will show how 
it imprints a signature of the KAM-Cantori transition and pro- 
vides experimental realization of various landmarks of non- 
linear systems such as Arnold tongues and bifurcations. 

A. Fragmented Fermi Sea 

Fermions in the extended phase have metallic properties. 
For A = 0, the single particle eigenstates are fully local- 



As disorder is increased, more quasi-Fermi seas become 
visible, and with the increase of v they overlap in a com- 
plicated pattern. Exactly at criticality, A = 1, the pattern 
evolves into a fractal-like structure. Beyond the critical point, 
A > 1, the fragmented quasi-momentum distribution profile 
disappears and instead it becomes a smooth function of Q and 
V. This behavior is summarized in Fig. |4]and the correspond- 
ing cross sections for fixed Q and v are displayed in Fig. |5] 
and Fig. |6] 

At irrational filling factors, however, the quasi-momentum 



ized in quasi-momentum space (and thus delocalized in po- 
sition space) and n((5) is a step-like profile: n(ff) = 1 for 
\Q\ ^ Qf andn((5) = for \Q\ > Qp (Qp is the Fermi mo- 
mentum). For < A < 1 the single particle eigenstates local- 
ized at k acquire some admixture of other quasi-momentum 
components ( nk±FM-i leading order in A). In this regime 
the quasi-momentum distribution retains part of the step-like 
profile but gets fragmented into additional structures centered 
at different reciprocal lattice vectors of the secondary lattice. 
We call the filled states centered around Q = the main Fermi 
sea and those around the QP related reciprocal lattice vectors 
the quasi-Fermi seas. 



The fragmentation of the quasi-momentum distribution and 
the development of QP Fermi seas can be understood from 
first order perturbation theory. To first order in A the momen- 
tum landscape becomes fragmented in six regions, shown in 
Fig [3] and given by: 



(12) 
(13) 
(14) 

(15) 

(16) 
(17) 

by white lines in 



(18) 



distribution profile remains smooth regardless of the value of 
A. The smoothness of the momentum distribution at these spe- 
cial filling factors can be understood using the same reasoning 
as the one used to understand the smooth character of the re- 
turn map in position space in terms of non-degenerate pertur- 
bation theory. In this case, however, the states that are coupled 
are the ones localized at the quasi-momentum Qf + eg and 
-Qf + eg, ( ega <C 1). 

The fragmentation of the Fermi sea in the A < 1 regime, 
and the smooth profile in both the A > 1 and band insulator 
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FIG. 3: (color online) Various boundaries of the fragmented Fermi 
sea. 



phases are directly connected, via the self-duality property, to 
the discrete nature of the density return map in the localized 
phase, and its smooth character in the extended and band in- 
sulator phases. 

An important point to emphasize again, which is crucial for 
possible experimental observation of the predicted behavior, 
is the insensitivity of the momentum distribution to variations 
of the phase (j). We demonstrate such insensitivity by noticing 
that the quasi-momentum distributions plotted in Fig |4] are ac- 
tually averaged over many different values of (j). Hereafter 
all the plots in momentum space are always averaged over 50 
random phases. 



B. Arnold Tongues 

Arnold tongues are mode-locked windows that characterize 
the periodic dynamics of iconic non-linear systems with com- 
peting periodicities. In the 17th century Christian Huyghens 
1 25 1 noted that two clocks hanging back to back on a wall 
tend to synchronize their motion. In general coupled systems 
such as coupled pendula or pendula whose lengths vary peri- 
odically with time exhibit mode-locking |26|. As the parame- 
ters of a system are varied, it passes through regimes that are 
mode-locked and regimes which are not. 

Arnold-like Tongues also appear in the QP fermionic mo- 
mentum distribution, where they reflect the complex nature of 
the physics induced by the competing periodicities. 

For A = 0, the Fermi distribution is a binary distribution 
with only two values n((5) = 0, 1. For finite A and for generic 
filling factors, two windows of values centered around zero 
and unity become allowed whose width increases with A as 
can be derived from a perturbative analysis. We call such win- 
dows QP tongues (See Fig. [7]). Exactly at the onset of the 
metal-insulator transition the two distributions overlap, mim- 
icking the Arnold tongue behavior. We have checked that the 
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FIG. 4: (color online) Momentum distribution for (a)A = 0.1, 
(b)A = 1, (c)A = 2. All are averaged over 50 random phases to 
mimic a realistic experimental case. We can see that the sharp fea- 
tures of these distributions survive phase averaging. 



duality of the Harper equation allows one to observe the for- 
mation of analogous structures in position space. The absence 
of a metal insulator transition at the irrational filling factors 
at which the system is a band insulator is also signaled by 
the disappearance of the Arnold tongues at these fillings, as 
shown in Fig.|7jb). 
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FIG. 5: (color online) Momentum distribution for fixed Q{k = 
1045) , Fm = 4181 as a function of filling factors vs. the disor- 
der parameter. 
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FIG. 6: (color online) Momentum distribution vs the disorder param- 
eter, (a) v = 0.5, and (b) v = Fm-i/ Fm- 



C. Bifurcations 

Bifurcations are common features observed in nonlinear 
dynamical systems, which occur when a small smooth change 
made to specific parameters (the bifurcation parameters) of 
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FIG. 7: Momentum distributions a function of disorder parameter, 
(a) filling factor v — 0.5. Here we see the formation of two Arnold 
tongues which merge at the critical value A = 1. (b) z/ = a. When 
V is irrational, the system is a band insulator and no Arnold tongues 
are present. We have checked the existence of the similar structure 
in position space when converting A to 1/A. 



the system causes a sudden 'qualitative' or topological change 
in its behavior. A series of bifurcations can lead the system 
from order to chaos. 

The density distribution provides a nice manifestation of a 
single bifurcation with A as a bifurcation parameter. When the 
density at the various Fibonacci sites is plotted as a function 
of A, for an specific filling factor which depends on the value 
of (j), a bifurcation opens up at A = 1. In Fig. [8] the existence 
of a bifurcation at quarter filling {v = 1/4) when (j) = 37r/4 
is shown. 

A qualitative understanding of the bifurcation can be gained 
by considering the two limiting regimes, A <C 1 and A 1. 
In the weak coupling limit, A <C 1, the local density is uniform 
and directly proportional to the filling factor u. On the other 
hand, in the strong coupling limit, a given site remains empty 
or occupied depending upon whether the onsite potential is 
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greater than or less than the Fermi energy. 



Fibonacci sites j = have similar on-site energies which oscillate about Cc = 2 A cos (^) (assuming n ^ M and (j) > 

27rFM-n/FM ) as: 



2Acos(27rcrFn + 0) = 2Acos (27r(-l)''"^FM-n/i^M 
These oscillations can be seen in FigfTSlin Appendix A. 



r 



< Cc n is odd 
> Cr n is even 



(19) 




FIG. 8: (color online) Bifurcation of local density at Fibonacci sites. 
The central curve (red triangle) corresponds to quarter filling while 
the lowest (black blocks) and the topmost (blue disks) curves re- 
spectively correspond to u = 1/5 and u = 1/3. Quarter fill- 
ing is a special case where the local density at Fibonacci sites is 
or 1 as A oo. For v > 0.25, the Fibonacci sites are filled 
while for < 0.25 Fibonacci sites are empty as A oo. We 
have checked the existence of the same bifurcation phenomena in 
the quasi-momentum distribution but with the weak and the strong 
coupling limits reversed ( A ^ 1/A). 

Consequently, at the filling factor, at which the Fermi 



energy matches Cc (e.g. if (j) = 37r/4 then Uc ~ 0.25) Eq. 19 
implies that as A goes to infinity, even Fibonacci sites (F21) 
become empty and odd Fibonacci sites (F2/+1) become oc- 
cupied. This behavior combined with the monotonic increase 
of the density with A in the weak coupling limit qualitatively 
explains the observed bifurcation. 

The self-dual behavior described by Eq[TT]implies the exis- 
tence of the same bifurcation phenomena in quasi-momentum 
distribution but with the weak and the strong coupling limits 
reversed ( A ^ 1/A). 

When the local density is plotted for a generic lattice site as 
a function of A and the landscape shows a similar change in 
topology. In contrast to Fig. [8] where the filling factor is fixed 
but different curves corresponding to different Fibonacci sites 



P377 




PlOO 



(a) 




(b) 



FIG. 9: (color online) Local density distribution as the disorder pa- 
rameter and filling factor are varied for (a): Fibonacci site j = 377 
and (b): normal site j = 100. We have checked the self-dual behav- 
ior in the quasi-momentum distribution with the weak and the strong 
coupling limits reversed ( A ^ 1/ A). 



are shown, in Fig. |9]the density at an specific lattice site is 
plotted but both the filling factor and disorder strength are var- 
ied. In this plot besides the large bifurcation there are smaller 
ones which cannot be explained by studying the two limiting 
cases. However, in momentum space they are qualitatively 
understood by considering cuts at a fixed A as a function of u. 
In these cuts the big jump at small A can be identified with the 
edges of the main Fermi sea in the quasi-momentum profile 
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and the smaller jumps correspond to edges of the quasi-Fermi 
seas. 



D. Devil's Staircases 

A Devil's staircase describes a self similar function /(x), 
with a hierarchy of jumps or steps. In other words, f{x) 
exhibits more and more steps as one views the function at 
smaller and smaller length scale. The derivative of f{x) van- 
ishes almost everywhere, meaning there exists a set of points 
of measure such that for all x outside it the derivative of 
f{x) exists and is zero. 

Cold atom experiments may provide an opportunity to vi- 
sualize devil's staircases in the momentum-momentum corre- 
lations, namely the noise correlations. 

Measurement of noise correlations is an example of 
Hanbury-Brown-Twiss interferometry (HBTI), which is sen- 
sitive to intrinsic quantum noise in intensity correlations. 
HBTI is emerging as one of the most important tools to pro- 
vide information beyond that offered by standard momentum 
distribution-based characterization of phase coherence. The 
noise correlation pattern in ID QP bosonic systems has been 
studied theoretically 1 27^-^301 and has also been measured ex- 
perimentally |22|. Here however we focus on the fermionic 
system. 

In the extended phase, noise correlations exhibit a series of 
plateaus as v is varied, and the number of steps or plateaus in- 
creases as the strength of the disorder A increases. The origin 
of this step-like structure with jumps occurring at the filling 



factors 



in) 



|(5n|<^/7r, can be understood from a perturba- 



tive argument as follows. For non-interacting fermions, for 

' ' 2 



Q ^ (See Fig, 



10),A(Q(fe)) = -E 



(m) (m) 

Tk To 



with 



being the Fourier transform of the m^^ single-particle 
eigenfunction. For A = 0, the overlap between any two dif- 
ferent Fourier components is always zero as only the ground 
(m = 1) state has a zero quasi-momentum component, i.e. 
1. For small A <C 1, first order perturbation the- 



ory yields a single step observed at u^]^ = |(5i|a/7r, as only 



with m*^^^ 



y^^} Ni is nonzero. 



Quantitatively, the heights of the steps at v^^^ for A <C 1 are 
given by: 



C(o) = 
C(Qi) = 



1 



A 

2) (1 - cos(27ra))2 
1 

(1 - cos(27rcr))2 



(20) 



Where C(0) (C(Qi)) is the step height for A(0) (A(Qi)) 
at filling factor v^^^ . The minus sign implies a decrease in 
the noise as the filling factor u increases. As A increases, 
more and more steps are seen and can be explained using 
higher order perturbation theory. At criticality (see FigsfTO 



11), the steps acquire a hierarchical structure which resem- 
bles a devil's staircase and which correlates with the fractal 



0,08 ^ 
0.06\ 

-A(Q)xl0^0'04 

0.00 



-n 
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(b) 



-71 




(c) 



FIG. 10: (color online) Noise correlation for (a) A = 0.1, (b)A = 1, 
(c)A = 2. All plots display averaged quantities over 50 random 
phases. We can see that structure in the noise correlations sur- 
vives phase averaging, as for the momentum distributions of Fig. 4. 
A(Q = 0) is not displayed in those plots. 



structure of the energy spectrum ||23l[24]|. In contrast to the 
momentum distribution, noise correlations do not show sig- 
nificant differences between the rational and irrational filling 
factors. 

In the localized phase, A > 1, we see the smoothing of the 
step- structure and the noise correlation function tends towards 

a sinusoidal profile A((5n) ^ ^(^n,o - '''\}n)'^^ • 
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To emphasize the striking difference between the noise pat- 
terns in the extended and locaHzed phases (the step-Hke vs 
smooth profile), in Fig. [12] we plot a measure of the gaps de- 



}.0 0.2 0.4 0.6 0.8 1.0 

y 

(a) 




'%.0 0.2 0.4 0.6 0.8 1.0 

y 

(b) 



FIG. 11: (color online) This figure shows the evolution of A{Qn) 
(n = 0, 1) as a function of the filling factor for different A values. 
In the limit A ^ 1, it exhibits steps occurring at the filling factors 
Uj^"^ (m — 1,2,...) (see text). In the limit A 1, A(Qn) acquires 
a sinusoidal profile. The inset shows the self-similar nature of the 
steps at A = 1. 



veloped in A(0) as a function of A. 



VI. SUMMARY 

We now summarize our key results: 
1 . Extended-localized vs KAM-Cantori Phases 

At the many-body level, the connection between the 
metal-insulator transition and the KAM-Cantori transi- 
tion is signaled in the return map of the local density 
between nearest neighbor sites. For a given angle 0, the 
metallic phase with extended single particle wave func- 
tions exhibits a smooth return map. We identify this 
with the preserved invariant KAM tori in phase space at 
weak perturbations in the corresponding classical sys- 
tem. The localized phase displays a discrete return map, 
which we identify with the remaining tori or Cantori 
outside the perturbative regime. This behavior is shown 
in Fig[T] In addition, for the particular filling factors at 
which the system becomes a band insulator, the return 
map remains smooth for any value of A (see Fig. [2]). 

Similar behavior can be observed in the corresponding 
momentum distribution return map. However, while the 
local density return map loses these characteristic fea- 
tures after averaging over (/), the momentum distribu- 
tion return map remains almost unaffected. The robust- 
ness of the momentum return map to phase variations, 
is ideal for the experimental visualization of the metal- 
insulator vs KAM-Cantori connection in cold atoms. 




FIG. 12: (color online) This figure shows a measure of the gaps in 
A(0) vs A for different system sizes, Fm- We deem that a gap occurs 
if |A(0)(A/'p + l)-A(0)(A/'p)| > e(wechoosee= 10"^^ here). We 
can see that with increasing M, the measure of the gaps approaches 
a step function with the step position at A = 1. 



2. Quasi-fractal structures and Arnold Tongues 

The introduction of weak disorder modifies the char- 
acteristic step-function Fermi- sea profile of the quasi- 
momentum distribution. Additional step-function struc- 
tures centered at different reciprocal lattice vectors of 
the QP lattice appear for A > 0. We refer to those 
structures as the "quasi-Fermi seas". With increasing 
A and the number and width of the various quasi- 
Fermi seas increase the fragmentation of the momen- 
tum distribution, turning it into a complex pattern. At 
A = 1 the fragmentation becomes maximal and the mo- 
mentum distribution evolves into a smooth profile as 
the system enters the localized phase (See figures [3][6]). 
The overlap of the various quasi Fermi seas as one ap- 
proaches criticality is reminiscent of the Arnold tongues 
overlap observed in non-linear systems, such as the cir- 
cle map, as they enter the chaotic regime L31J . 

A more appropriate analogy of such behavior can be ob- 
served by considering the set of values taken by quasi- 
momentum distribution for a given filling factor. In 
the absence of disorder this distribution can only take 



11 



the values 1 or 0, depending upon whether the quasi- 
momentum is greater or smaller than the Fermi quasi- 
momentum. As the strength of the quasi-periodic lattice 
increases, two distributions of values develop, centered 
around and 1 respectively. Their width increases with 
increasing disorder and they overlap exactly at critical- 
ity (Figig. 

3. Bifurcations: The overlap between Arnold tongues at 
Ac, can be linked, using the space-momentum duality 
transformation, to the appearance of a bifurcation in the 
density profile. The bifurcation occurs at a common 
but phase dependent filling for the various Fibonacci 
sites (Fig[8|. At generic lattice sites, the filling factor 
at which the bifurcation takes place also depends on 
the lattice site under consideration and can be observed 
when the local density is plotted as a function of A and 

y (Fig|9l). 

4. Devil's staircases: Noise correlations plotted as a func- 
tion of the filling factor exhibit step-like structures 
which evolve into a devil's staircase at the onset to the 
metal-insulator transition ( figures 70p2 ). 



we use to link the quasimomentum-position observables. In 
the plot we highlight the Fibonacci sites. A small system size 
is used to make the visualization clearer. 



Systems with competing periodicities stand in between pe- 
riodic and random systems. The richness and complexity un- 
derlying such systems have been studied extensively ifTSl . Ul- 
tracold atoms are emerging as a promising candidate to sim- 
ulate a wide variety of physical phenomena. Here we have 
shown they offer opportunities to experimentally realize vari- 
ous paradigms of nonlinear dynamics. 

Our focus here was on spin-polarized fermionic systems, 
since we wanted to look at the simplest consequences of 
many -body physics in disordered systems. However, Bose- 
Einstein condensed systems may also be used as tools for 
laboratory investigation of various predictions made for the 
quasi-periodic systems based on single-particle arguments 
|[T6l . For example, it might be possible to confirm the strong 
coupling universality prediction, which establishes that the ra- 
tio of the single particle density at two consecutive Fibonacci 
sites should be a universal number |32|. 

Acknowledgments A. M. Rey and S. Li acknowledge sup- 
port from the NSF-PFC grant and NIST. 



Appendix A: Mapping from position to momentum space from 
self-duality relationship 

Here we provide an example of the mapping between posi- 
tion coordinates j and quasi-momentum coordinates k which 



80 
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£^40 
20 
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7 

FIG. 13: (color online) Mapping from j to k used to determine 
momentum- space observables from position- space observables and 
vice versa. Larger points are Fibonacci sites. Here we choose 
Fm — 87, and is the relation that maps j to k. 

Appendix B: Perturbation Theory and Dimerized States 

We begin our analysis with the Harper equation (|2]). For 
A oo, the single particle wave functions are local- 
ized at individual lattice sites, V^j^^ = (^j,L^ and e*^^^ = 
cos(27r crL^ + 0), where Lm is defined by cos(27rcrI/^_i + 
0) < cos(27rcrL^ + 0). When A ^ 1, we can get the sin- 
gle particle wave function through exact Diagonalization, or 
we can use the perturbation theory to obtain the approximate 
eigenvalues and eigenf unctions. 

We focus first on the site L^. Assuming that 
2A| cos(27rcrI/^ + (/))— cos(27rcr(L^ ± 1) + (/))| ^ 1, non- 
degenerate perturbation theory can be used, and the results 
are: 
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,(m,0) 



2A cos(27r crL„ 


1 



2A \cos(27rcrL^ 

1 



2A \cos(27r<jL^ 



1 



1 



- cos(27rcr(L^ + 1) + 0) cos(27rcrL^ + 0) — cos(27rcr(L^ — 1) + 0) 



cos(27r<j(L^ + 1) 



cos(27r(jL^ 



cos(27r<j(L^ — 1) 



(Bl) 
(B2) 

(B3) 
(B4) 
(B5) 



where e*^^'^\ e^^^'^^ ^(^'2) ^j,^ ^j^^ zeroth, first and second 
order terms of the eigenvalue e^^^ , with a similar meaning for 
andV^f Also, = ^l^,, + 0((^)2). So, up 

to first order of ^, J^^. IV^]"^^ p = 1. 

For other sites that satisfy 1 2A| cos(27rcrL^ + — 



cos(27rcr(L^ + 1) + (/))|, non-degenerate perturbation theory 
doesn't work and degenerate perturbation theory will be used. 
In that case, the zero order energies e and eigenfunctions 
satisfy: 



2Acos(2wi:. + </.) -1 U 4"°! V e f ^ (B6) 

The above equations result in a pair of energies which we denote as and Cd' 

Cd = X{cos{27TaLm + 0) + cos(27ra(L^ + 1) + 0)} - vTTA^ (B7) 

= A{cos(27raL^ + 0) + cos(27rcr(L^ + 1) + 0)} + VTTA4 (B8) 



where = A(cos(27r(jL^ + 0) — cos(27r<j(L^ + 1) + 
The corresponding orthogonal eigenfunctions can be written 
as: 



and = 



B 



-A 



(B9) 



where A 
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. + 1 ' V ^v'A^ + i 
choose a (/3) to be the greater (lesser) of A and B\ \a 

yi+A^ + IA^I \o\2 Vl+A^-|A.^| 



We can 



2V1+A4 



2V1+A2 



From these results one 



can see that \af increases as |A^| increases. For example, 
when |A^| = 0, |ap = 0.5 and when |A^| = oo, |ap = 1. 

When |A^| <C 1, is comparable to |ap and the states 
are localized at the same two neighboring sites. Hence, they 
will be referred as a pair of dimerized states ( The numerical 



results are shown in Fig 14 (a) ). The dimerized states are 
found by noticing that: 



cos(27rcrL^ -\ 
2 sin {27T(j{Lr, 



I — cos(27rcr(L^ + 1) 
1/2) + 0) sin(7rcr) 



(BIO) 



So for the points that satisfy | sin (27r cr(L^ + ^) + (/)) | <C 
Y, dimerized sites exist at Lm and Lm + 1. Since a = , 



it is clear that we can construct the relation Mod[27rcrL^ + 
0,27r] ^ 27r Mod[5^/FM,l], where m G (1,2,...,Fm) 
are indexes of increasing energy level, Lm G (1,2,..., Fm) 
are indexes of sites' positions and Bm = Mod[L^ x Fm-i + 
^ Fm] ^ (0, 1, 2, . . . , Fm — 1) are indexes that we introduce 
for the convenience of discussion (/ is the closest integer to 
^Fm)- Using this relation it can be shown that m ^ 2\Bm — 
Fm/2\ (we use ^ here, since the exact expression for even 
and odd Fm are slightly different). The relations of m, Lm 
and Bm are plotted for the case of Fm = 610 in Fig [14] (b). 
As we discussed, the dimerized states Lm satisfy: 



sin(27rcr(L^ + 1/2) + 0) 



Writing the above equation in terms of Bm, we get: 



2lTBm 1 ^ 

sm I —- h TTcr I ^0. 

^M 



The solutions are ^^"^^.^^ ^ ^ 1 or 2. Therefore, 



(Bll) 



Fm 



Br, 



^mi ~ {Fm — Fm-i)/2 

Bm2 ^ (2i^M — Fm-i)/2 



(B12) 
(B13) 
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(b) 



FIG. 14: (Color online): Panel (a) shows numerically that at the two 
red points, mi = Fm-i and 1712 = Fm-2, Am vanishes, so that 
degenerate perturbation theory is required. Panel (b) shows the rela- 
tions of Lm, Bm and m (see text) for = ^ . The relations between 
Lm and m are very sensitive to the value of (j) while the relations be- 
tween Byn and m are not sensitive at all. 



which correspond to: 



mi ^ \2{B^, - Fm/2)\ - Fm-1 (B14) 
m2 ^ |2(5^, - Fm/2)\ ^ Fm-2 (B15) 



A detailed analysis shows that the band opens exactly at 
mi = Fm-1 and m2 = Fm-2 regardless of Fm being even 
or odd. This demonstrates, at perturbative level, the special 
behavior of the return map at irrational filling number Fm-2 
orFM-i(See Fig^. 

In our calculations, we find out those Lm with the low- 
est values of | sin (27rcr(I/^ + 1/2) + 0) | (avoiding double 
counting of in different pairs) and then apply to those 
points the perturbation theory we discussed above. 

Based on the latter considerations it is possible to demon- 
strate that around mi = Fm-i or 1712 = Fm-2, there exists 
a sequence of paired states, exhibiting the following relation- 
ships as A ^ 00, 



with i <C Fm and r = 1,2. These states determine the prop- 
erties of the system when the Fermi energy is close to the ma- 
jor gaps, i.e. u = a or u = a'^, and thus the main properties 
of the band insulator phases. 
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